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Ripples and Ripples: from Sandy Deserts to Ion-Sputtered Surfaces 
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We study the morphological evolution of surfaces during ion sputtering and we compare their 
dynamical roughening with aeolian ripple formation in sandy deserts. We show that the two phe- 
nomena can be described within the same theoretical framework. This approach explains the dif- 
ferent dynamical behaviors experimentally observed in metals or in semiconductors and amorphous 
systems. In the case of ion erosion, we find exponential growth at constant wavelength up to a 
critical roughness Wc- Whereas, in metals, by introducing the contribution of the Erlich-Schwoebel 
barrier, we find a transition from an exponential growth to a power law evolution. 

PACS numbers: 45.70.-n 61.43.-j 68.35. Ct 68.55.Jk 



I. INTRODUCTION 

When an ion hits a surface liberates locally a large 
amount of energy that melts a region of the solid im- 
mediately below. For geometrical reasons, the sputter- 
ing effect depends on the surface-curvature: the energy 
concentrates on regions of positive curvature and this fa- 
vorites the excavation of valleys and the growth of hills. 
On the other hand, thermal diffusion and surface ten- 
sion tend to smoothen the irregularities by flattening the 
surface. It has been observed that under the combined 
action of these mechanisms the surface tends to create 
spontaneously ripples In nature, 

ripples are commonly observed in sandy-deserts as the 
result of dynamical instability of the sand surface under 
the action of a sufhciently strong wind In this case, 
the formation of ripples is commonly associated with the 
effect produced by some grains that are lifted from the 
sand-bed and accelerated by the wind. These grains, 
when re-impact with the bed, splash up a number of other 
grains. Most of these grains return to the bed leading to 
a local rearrangement, whereas some other are acceler- 
ated by the wind and impact again on the bed after a 
certain 'saltation' length. In the literature, many studies 
have been devoted to understan ding the mechanism of 
ripple formation [13 [H El El El El- In particular, 
an hydrodynamical model for aeolian ripple formation, 
based on a continuum dynamical description with two 
species of grains (immobile and rolling grains), was pro- 
posed with success by Bouchaud et al. [S El El Hi- 
The main ingredient of such a model is a bilinear differ- 
ential equation, for the population of the two species of 
grains, which shows the instability of a flat bed against 
ripple formation. 

In this paper we show that the same reasoning which 
has been used to describe the sand ripples formation 
in deserts, applied to the studies of dynamical surface 
roughening, leads to an accurate description of the mor- 




FIG. 1: Ripples on sand (Gobi desert) and Ripples on sur- 
faces (Ag under ion sputtering). 



phogenesis and evolution of ripples on crystal and amor- 
phous surfaces during ion sputtering. The present ap- 
proach contains not only the Bradley-Harper approach 
|20l I2II (based on a Kardar-Parisi-Zhang type equa- 
tion [231 ), but it is also able to describe some of the 
crucial experimental features recently observed in these 
systems In particular, by means of this approach 
we can explain the two distinct dynamical behaviors 
experimentally observed in amorphous/semiconductors 
systems and in metals 1^. In the first case (amor- 
phous/semiconductors) we find that the ripples growth 
exponentially fast at constant wavelength A up to a crit- 
ical roughening Wc at which the growing process inter- 
rupts. On the other hand, in metals (when the Erlich- 
Schwoebel barrier is active) , we find a transition between 
an initial exponential to a slower power-law growth of the 
roughness (the root mean square of the height profile). 
In this regime the ripple- wavelength tends also to growth 
with time. 



II. PARTICLE MOBILITY AND RIPPLE 
DYNAMICS 
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When the surface of a solid is taken under ion sput- 
tering some atoms in the proximity of the surface receive 
energy from the sputtered ions and pass from a bounded - 
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'immobile' - solid state to a 'mobile' melted unbounded 
state. The opposite mechanism is also allowed: some 
mobile atoms can gain in energy by becoming immobile 
and bounding in a given position in the solid. A certain 
fraction of atoms might also be dispersed into the atmo- 
sphere. Let us call h{r, t) the height of surface profile 
made of immobile -bounded- atoms and call i?(r, t) the 
height of mobile -melted- atoms. In analogy with the 
theory developed to explain the dynamical evolutions of 
dunes in deserts 0, ' describe the mech- 

anisms of excavation, exchange between mobile and im- 
mobile atoms and surface displacement of mobile atoms 
in term of the following differential equation: 

dh 

— = -r{B,h),,+r{R,h)ad 
at 

3R 

— = ~yJiR,h) + {I - ^)T{R,h),, ~r{R,h)ad (.1) 

Where T{R, h)ex is the rate of atoms that are excavated 
under the action of the sputtering, and (1 — </)) is the part 
of them that pass from immobile to mobile, whereas (f) is 
the fraction that is dispersed into the atmosphere. 

Let us now write in details the various terms contained 
in EqH 

A. Excavation 

The excavation effect must clearly depend on the num- 
ber and velocity of the sputtered ions (i.e. its flux), but 
also the local shape and orientation of the surface might 
play an important role. Indeed, the energy transmitted 
by the impacting ions concentrate more in regions of the 
surface with positive curvature. Moreover, part of the 
surface facing the flux are likely to experience a different 
erosion respect to others which are less exposed to the 
flux. Crystalline orientation and anisotropics might be 
also taken into account. We can write: 

r(i?,/i)e:r = +aV/i + 5V2/i) ; (2) 

here 77 is the sputtering flux, whereas a and b are respec- 
tively associated with the flux-direction-dependent and 
with the curvature-dependent sputtering erosion. 

B. Adsorption 

The rate of adsorption of mobile atoms into immobile 
solid positions must be dependent on the quantity of mo- 
bile atoms in a given spatial position. Similarly to the 
excavation process, the adsorption is also dependent on 
the local curvature and orientation. We can write: 

T{R,h)ad^ Ril + cVh + dV^h) , (3) 

where the parameter 7 is the recombination rate and c 
and d are associated to the different probabilities of re- 
combination in relation with the local orientation and 
shape of the surface. 



Note that Eqs. |51 U] contains the same terms as 
the ones proposed in the literature for the formation 
of aeolian dunes in the so-called hydrodynamical model 
[Til Ir^ ITM 0, |2^. Indeed, in deserts, sand grains are 
lifted from the sand-bed and readsorbed into it with a 
probability which is dependent on the local shape and 
orientation of the dunes. Eqs. [3 El represent the sim- 
plest analytical expressions which formally take into ac- 
count these shape and orientation dependences. In the 
search for simple explanations, such equations are there- 
fore rather universal. 



C. Mobility 

Mobile atoms will move on the surface, and the quan- 
tity 3{R,h) in Eq^is the 'current' of these atoms. In 
surface growth, there are two main mechanisms that are 
commonly indicated as the responsible for the surface 
mobility of atoms The first is a current, driven 

by the variations of the local chemical potential, which 
tends to smoothen the surface asperity moving atoms 
from hills to valleys. The second is the current induced 
by the Erlich-Schwoebel barrier which -on the contrary- 
moves atoms uphill. In addiction to these main mech- 
anisms we might also have to take into account a drift 
velocity and a random thermal diffusion, obtaining: 

3iR,h) = lCRS/{V^h) + sR ^^^——r + wR-VVR . 

1 + (adyhy 

. . 

In this equation, the first term describes a deterministic 
diffusion driven by the variations of the chemical poten- 
tial which depends on the local shape of the surface; the 
second term is associated with the uphill current due to 
the Erlich-Schwoebel barrier and ad is a constant associ- 
ated with the characteristic length. The quantity v is a 
drift velocity of the mobile atoms on the surface, whereas 
V is the dispersion constant associated with the random 
thermal motion. 

Note that Eq0|is substantially different from the one 
proposed in the literature to describe ripples in granular 
media [IE 113, EE ES 11^ • Here the current is supposed 
to be dependent on the local shape and orientation of 
the surface (the h{r,t) profile). The equations describ- 
ing sand deserts can be retrieved from Eq^ by impos- 
ing IC — and s = 0, but -on the contrary- in surface 
growth these two parameters are the leading terms of the 
equation and play the role of control parameters in the 
dynamics of ripple formation. Nonetheless, these terms 
describe a rather simple dependence of the dynamic of 
particles on a surface on the geometrical shape of the 
surface itself. Again, in our seek for universality, we ex- 
pect that similar terms can be profitably introduced in 
the context of aeolian sand ripples in order to describe 
specific phenomena (associated, for instance, with pack- 
ing properties or granular flow pE|) which relate the 
current of grains with the dune-shapes. 
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FIG. 2: Numerical solutions of EqQat various times indicate 
that under the action of ion sputtering the surface develops an 
instability which leads to the formation of ripples with a well 
defined characteristic wavelength. In the figure the black-tick 
line is the final surface-profile, whereas the tinnier-gray lines 
are some profiles at previous times. 



It should be noted that the factors a, c and v in Eas l2l3l 
and 21 are vectors (i.e. they have -in general- different 
components in the two horizontal directions). Indeed, 
crystal surfaces are in general anisotropic and therefore 
one must take into account the dependence of the param- 
eters on the relative orientation of the crystal-surface and 
the sputtering direction. 



III. DISPERSION RELATION 

A trivial solution of EqQ can be written for a com- 
pletely flat surface: h{r,t) = ho{t) and R{r,t) = Rq. 
In this case, we obtain Rq — (I ~ (t>)ri/j and ho{t) = 
—(j)rit+const.. This describes a surface that rests flat and 
it is eroded with a speed equal to (j)?]. But this behavior is 
only hypothetical since -in general- the dynamics of the 
surface-profile presents instabilities against spontaneous 
roughening and therefore its evolution is more complex. 
For instance, a numerical solution of Eq^ is shown in 
Fig|21(for the 1-dimensional case). We observe that, in a 
certain range of the parameters, the surface is unstable 
and periodic ripples are formed spontaneously. 



A. Stability analysis 

In order to infer indications about the amplification or 
the smoothing of small perturbations and to deduce an 
analytical expression for the ripples wave-length at their 
beginning, we performe a stability analysis on Eq^ For 
this purpose we assume that the surface-profile is made 
by the combination of a flat term plus a rough part: 



R{v,t) 

h{x, t) 



Ra + Ri{r,t) 
ho{t) + hi{r,t) 



with i?i(r, i) = Riexp{iujt + zkr) and hi{r,t) = 
/ii exp(iwt -|- ikr). We substitute these quantities into 
Eq^ and linearize the equation by taking only the first 
order in i?i and hi. A Fourier analysis (see Appendix 
shows that such a linearized equation admits solutions 
when the frequencies uj and the wave vectors k satisfy: 

[iuj + 'j + ikv + k^V] ■ 

{lOJ + ik [Vi - (1 - ,^)V2] - [Vi - (1 - 4>)V2] ] - 

7(1 - 0)[ik(vi-V2)-A;2 {V1-V2 - .si)-fc*/Ci] = ; 

(6) 

where, to simplify the equations, we have introduced the 
following notation: 

vi = 77a V2 = 770/7 

Vi =r]b P2 = r]d/j 

si = 7]s/j tCi = ?7/C/7 

Equation El establishes a dispersion relation w(k) that is 
a complex function with two branches corresponding to 
the solutions of the quadratic Equation (HI 



IV. SURFACE INSTABILITIES 

The kinetic growth of the surface instability is related 
to the immaginary part of Lo{k). Indeed, 7m(a;(k)) corre- 
sponds to modes with amplitudes that change exponen- 
tially fast in time, and negative values correspond to un- 
stable modes that increase with the time. We can there- 
fore study Imiuj) from the solution of EqjSl and search 
for the range of k in which Im{u}) is negative. The most 
unstable mode is the one that grows faster and it corre- 
sponds to the value of k at which Im{uj) reaches its most 
negative value (see Figj^l). 

The solution of EqjSlfor Im{uj), is 



2Imiuj)± 7+ V-Vi + {1- (j))V2 



± 



Ai-f (A2-f 4A2)i/2 



(7) 



where we have 



Ai = 7^ - { V - vi + (1 - 0)v2 k| 



27 



P - (1 - 20)Pi + (1 - (f>)V2 + 2(1 - 0)si 



and 



A, 



+ { [(P + Pi - (1 - 0)^2] ' - 47(1 - 0)/Ci } 



7 V + (1 - 20)vi - (1 - (f>)v2 k 

P + Pi - (1 - (/))P2j [v - Vi + (1 - 0)V2 



(8) 



(5) 



(9) 
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FIG. 3: The imaginary part of the dispersion relation Im{u)) 
can assume negative values which are associated with the sur- 
face instability (arbitrary units) |29||. The amplitude of modes 
with wavelengths A > 27r/fc* will grow exponentially fast. The 
tick line is the imaginary part of the analytical solution of 
Eq|H] whereas the tinny-gray line is the approximated expres- 
sion (at the fourth order in k) obtained for small ion flux (77 
small) . 



Let us first observe that in absence of sputtering (i.e. 
when ry = and therefore, Vi = 0, V2 = 0, Vi = 0, 
I?2 = 0, si = 0, ICi = 0) the solutions of Eq|Sl are 
w(k) — and w(k) = — kv + 1(7 + k^V). In this case, 
the imaginary part of w(k) is non- negative, therefore we 
-correctly- expect no spontaneous corrugation of the sur- 
face. On the contrary, when the sputtering is active 
(ry ^ 0), the immaginary part of a;(k) can assume nega- 
tive values as shown in Fig|31 where a plot of /m(w)_ is 
reported (along one direction of the vector k) . As one can 
see, typically the branch Im{Lo)- takes negative values 
for k between and a critical value k* at which it passes 
the zero |29| . The critical point k* , fixes the minimal un- 
stable wavelength. We therefore expect to find unstable 
solutions associated with the formation and evolution of 
ripples with wavelengths X> \* =2n/k* . 



V. RIPPLE WAVELENGTH 

Several analytical solutions of Eq^l can be found in 
special cases which are discussed in Appendix^ But the 
study of the surface instabilities can be highly simplified 
if we consider the first order effects when the sputtering 
fiux rj is small. 



A. Approximate equation 

In the case of small spattering fluxes, the branch of 
/m(ti;(k)), with negative values can be approximated to: 



Im{uj). 
with 



Plfc6 + P^k^ + P3(k)fc2 + P4fc2 + Pj.(k) 



D2fc4 + 27X>fc2-f(vk) 24-72 



(10) 



Pi = V[{1 - 0)7/Ci - P(Pi - (1 - 0)P2)] 

P2 = (1 - 4>)-i[viV2 - si) + 7/C1] - (1 + 4>)-iVVi 

P3(k) = -[Pi + (l-<^)p2](vk)2 

P5(k) = (1 - 0)7(vk)[(vi - V2)k] (11) 

When k is sufficiently small ( A; ^ 7/77 ), we can de- 
velop EqEIIat the 4*'* order obtaining: 



/m(w)_ ~ Ak'^ - Bk^ , 



(12) 



with 



4- ^2 



1 



, . 2jV + u2 
v{vi - V2) ^ 



B = (l)Vi + {l 



Sl 



v{vi - V2) 



(13) 



Here v, vi and V2 are respectivelly the components of v, 
vi and V2 in the direction parallel to k). (In FigO a 
comparison between this approximate solution and the 
exact one is given.) 



B. Solutions 

The expected wavelength of the ripples is associated 
with the fastest growing mode, which corresponds to the 
value of k at which Im(uj)- reaches its most negative 
point. Here the minimum of Imiuo) is at 




(14) 



Therefore, at the beginning, the roughness will grow ex- 
ponentially fast asW ^ exp(p2i/(4^)) with associated 
ripple- wavelength at: 




(15) 



C. Some special cases 

Let us first observe that, when Ki, si and (j) — Q, the 
ripple wavelength, given by Ea ll5l coincides with the one 
found for sand dunes in deserts (see for instance [l^). 
In our notation the 'reptation length' is Iq = w/7, the 
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FIG. 4: Evolution of tlie surface rougiiness (log-scale) v.s. 
time (linear scale) from numerical simulations (see Appendix 
O . (Arbitrary units). The insertion is the log-log plot of the 
last part of the evolution. When the Erlich-Schwoebel barrier 
is active the dynamical evolution of ripples is characterized 
by an exponentially-fast growth at the beginning and then 
a 'saturation' to a power-law growth (linear trend in log-log 
scale) . 



'cut-off length' is Ic = (1^2 — T^i)/v, whereas vi — V2 is 
the collective drift velocity of the dunes. The approxi- 
mations usually applied in this context [T^ H^. imply: 
Ic ^ y/vjj, and jlc ^ f 1 — V2. Giving, from EatTsI 



A - 27r 



2vIqIc 



(16) 



Let us now consider the dynamical evolution of a 
surface under ion sputtering and in particular the case 
when the effect of the Erlich-Schwoebel barrier is not 
present (as for semiconductors and glasses). In this case, 
s = 0, si = and we also expect that the drift veloc- 
ity V and the dispersion constant V are equal to zero or 
infinitesimally small. Indeed, here the current of mobile 
atoms on the surface is mainly induced by the differences 
in the chemical potential. Under these assumptions, from 
EallSI the wavelength of the most unstable ripple is: 



(17) 



where we called v = ^h<f>/ {1 — 0), a quantity which plays 
the role of an effective surface tension. Note that Equa- 
tionElis the same result as from the Bradley and Harper 
theory US il IH US HH . 

When the Erlich-Schwoebel barriers are active (s, si ^ 
0), effects can be observed on the ripple- wavelength at 
their beginning, which becomes: 



A - 27r 



2/C 



VI. EXPONENTIAL/POWER-LAW GROWTH 
AND CRITICAL ROUGHNESS 



In metals, when the Erlich-Schwoebel barrier is ac- 
tive, there is an important non-linear contribution in 
the current of mobile atoms which becomes sizable when 
the roughness becomes sufficiently large and therefore 
((adV/i)'^) '-^ 1 (see EqEJ where the average is over all 
the surface positions. We observe numerically that this 
changes the law of growth of ripples: from an exponential 
to a power-law behavior. This effect is shown in Fig^ 
Numerically, all the computed exponents follow in the 
range between 0.65 and 0.85. The theoretical evalua- 
tion of this exponent is under current investigation. In 
this regime the ripple- wavelength tends also to grow with 
time. Experimentally, power law growth of the roughness 
and growth of characteristic wavelengths were observed 
in erosive sputtering on Ag(OOl) 4;j. 

In semiconductors or glasses, when no Erlich- 
Schwoebel barrier is present, it is physically intuitive that 
the exponential growth of the surface roughness (which 
is a characteristic of the beginning of the surface insta- 
bility) cannot continue indefinitely. Indeed, from the 
expression R{r,t) = Rq + Ri exp{iujt + ikr), which we 
used to derive Ea ll5l we can immediately observe that 
when ^1 > i?o = (1 — (j))ri/"f, the amount of mobile 
atoms might become negative. Since a negative amount 
of atoms is physically impossible, the process of exponen- 
tial roughness-growth described above must necessarily 
finish around a critical roughness given by: 



1 



(19) 



(18) 



This behavior is confirmed by numerical solutions of Eq^ 
and it is expected to be observable in semiconductors and 
glasses after sufficiently long times. 



VII. CONCLUSIONS 

We have shown that the same theoretical approach in- 
troduced to describe the formation of aeolian sand ripples 
can be conveniently applied to the study of the formation 
of periodic structures on surfaces under ion sputterning. 
Although the two phenomena are rather different, they 
can be described within the same conceptual framework 
by using rather general ideas that relate mobility, exca- 
vation and adsorption rates with the surface shape and 
orientation. 

We have obtained general expressions for the ripples 
wavelength in term of the system parameters. It has 
been shown that in some particular cases such a solution 
coincide with the ones already known in the literature for 
sand dunes and surface instabihty [T7LITll20ll2lll2^. We 
have discussed the effect of the Erlich-Schwoebel barriers 
and compared the result with numerical solutions. We 
pointed out the Erlich-Schwoebel barrier can be respon- 
sible for a dramatic change in the system dynamics: from 
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the exponential growth to a power law. Finally, we have 
demonstrated that the occurrence of a critical roughness 
is predictable within the present theoretical framework. 

It should be noted that the main purpose of this paper 
is to point out a relevant example of universality: two 
processes which have completely different scales present 
a dynamical evolution which obeys to the same geomet- 
rical constraints and thus can be described by using the 
same phenomenological model. On the other hand, we 
must observe that the class of solutions of Eq^is rich and 
complex - even in the linear approximation. Exhaustive, 
systematic studies of the classes of solutions of this equa- 
tion and their dependence on the set of parameters will 
be the subject of future studies and publications. 



APPENDIX A: FOURIER TRANSFORM OF THE 
LINEARIZED EQUATION 



where v, vi and V2 are the components of v, Vi and V2 
in the direction of k* . 

On the other hand when, /Ci, Si, I?i and P2 = 0, we 
find 



k* 



7(1; - ^Wi) 



V{vi - V - {1 ^ (j))v2) 



(B2) 



The effect of the deterministic diffusion induced by the 
chemical potential can be studied from the solution 



k* 



(1 - 0)7/Ci - V{Vi - (1 - cj))V2) ' 



(B3) 



which holds when v = 0, vi = 0, V2 = 0, si = and 
V-Vi + {1- (j))V2 > 0. 

Whereas, when vi, V2, "Di and 152 = 0, we find 



By substituting EqsEl H 01 and H into EqU and by 
neglecting the second order terms (in Ri and hi), we 
obtain the following linearized equation: 



dhi 
'dt 

dRi 



jRi - [vi - (1 - 0)v2] Vhi 
[Pi - (1 - 0)P2] V^/li 
-7i?i - vVi?i + PV^i?! + 

(vi - V2)Vft-l + 



(1- 

(Pi - V2 

A Fourier analysis of Eg I All leads to 
7^1 - |iw -I- ik [vi - (1 - 0)v2] - 
A:^ [Pi - (1 - ,^)P2] = 
[iuj + j + ikv + k'^V] 7?i - (1 - 0) 



(Al) 



ik(vi - V2) 



k"^ (Pi - P2 - si) - k'^K.i 



(A2) 



with i?i and hi the Fourier components of Ri and hi 
respectively. This equation is a simple linear equation 
in two variables. It admits a non-trivial solution when 
the determinant of the coefficients is equal to zero. This 
leads to EqEJ 



APPENDIX B: EXACT SOLUTIONS 

Analytical expressions for the value of k at which 
Im{uj) = (k*) can be calculated from Eq[7| in some 
special cases. 

In particular, when = 0, si = 0, /Ci = and P = 0, 
we obtain 



7(t;i - 1^2) 




{v - Vi + V2){V2 — Pi) 



(Bl) 



(B4) 



which implies that the uphill current due to the Erlich- 
Schwoebel barrier can generate instability even when the 
shape-dependent erosion and recombination terms are in- 
active. 



APPENDIX C: NUMERICAL SOLUTIONS 

The numerical solutions of Eq^ presented in this pa- 
per and in particular the ones shown in Figs|21 and 
21 have been performed as follows. We considered a 
one-dimensional flat substrate (/i(x, 0) = hO) of length 
L, with periodic boundary conditions. An infinites- 
imal quantity of mobile atoms were added randomly 
to the substrate (with < r(x, 0) < LIQ-^"). We 
then computed the profile-evolution using Eq^ with 
the derivative substituted with finite differences. To 
this purpose, substrate has been divided into N dis- 
crete points. The -adimensional- time indicated in Fig0| 
is the number of numerical steps. The height is in 
unit of L/N and the roughness is defined as w{t,L) = 

{[h{x',t) - {h{x,t))J'^)li^ (see, for instance, [H). 

Several computations with a number of points equal to 
TV = 100, 200 and 300 (the one pubfishcd have N = 300) 
have been performed to verify the effect of boundary and 
discretization. Moreover, simulations with no periodic 
boundary conditions and with the sputtering term (EqEJ 
applied only to a central mask, have also been performed 
obtaining very similar results. The robustness of the 
present approach has been verified varying the parame- 
ters, the time steps, the initial roughness of the substrate, 
etc. Comparable results have been always found but, we 
must stress that, under some conditions, numerical in- 
stabilities (in particular small surface-deformations with 
A ~ L/N) can be trigged on depending on the protocol 
utilized. 
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